#include <stdio.h>
#include <stdlib.h>
#include <complex.h>
#include <fftw3.h>

#include <math.h>

void fft_r2r(double *in,double *out,int size)
{
    fftw_plan_with_nthreads(4);
    fftw_r2r_kind kind = FFTW_REDFT00; // DCT-1
    fftw_plan plan = fftw_plan_r2r_1d(size,in,out,kind,FFTW_ESTIMATE);
    if(plan == NULL) {
        printf("fft:Create fftw_plan error!");
        exit(1);
    }
    fftw_execute(plan);
    fftw_destroy_plan(plan);
    fftw_cleanup_threads();
}

void fft_r2c(double *in,fftw_complex *out,int size)
{
    fftw_plan_with_nthreads(4);
    fftw_plan plan = fftw_plan_dft_r2c_1d(size,in,out,FFTW_ESTIMATE);
    if(plan == NULL) {
        printf("fft:Create fftw_plan error!");
        exit(1);
    }
    fftw_execute(plan);
    fftw_destroy_plan(plan);
    fftw_cleanup_threads();
}

void fft(fftw_complex *in,fftw_complex *out,long size)
{
    fftw_plan_with_nthreads(4);
    fftw_plan plan = fftw_plan_dft_1d(size,in,out,FFTW_FORWARD,FFTW_ESTIMATE);
    if(plan == NULL) {
        printf("fft:Create fftw_plan error!");
        exit(1);
    }
    fftw_execute(plan);
    fftw_destroy_plan(plan);
    fftw_cleanup_threads();
}

void ifft(fftw_complex *in,double *out,long size)
{
    fftw_plan_with_nthreads(4);
    fftw_plan plan = fftw_plan_dft_c2r_1d(size,in,out,FFTW_ESTIMATE);
    if(plan == NULL) {
        printf("ifft:Create fftw_plan error!");
        exit(1);
    }
    fftw_execute(plan);
    fftw_destroy_plan(plan);
    fftw_cleanup_threads();
}

void write_data(const char *filename,double *data,long size)
{
    FILE *pf = fopen(filename,"w");
    if(pf == NULL){
        perror("write_data:fopen()");
        exit(1);
    }

    for(int i = 0;i < size;++i){
        fprintf(pf,"%f\n",data[i]);
    }

    printf("Write successfully!\n");
    fclose(pf);
}

void write_data_c(const char *filename,fftw_complex *data,long size)
{
    FILE *pf = fopen(filename,"w");
    if(pf == NULL){
        perror("write_data:fopen()");
        exit(1);
    }

    for(int i = 0;i < size;++i){
        fprintf(pf,"%f\n",sqrt(creal(data[i])*creal(data[i]) + cimag(data[i])*cimag(data[i])));
    }

    printf("Write successfully!\n");
    fclose(pf);
}

void shift_freq(fftw_complex *freq,long size)
{
    fftw_complex temp;
    for(int i = 0;i < size / 2;++i){
        temp = freq[i];
        freq[i] = freq[i+size/2];
        freq[i+size/2] = temp;
    }
}

int main()
{
    int n = 64;
    // double * in = malloc(sizeof(double) * n);
    // fftw_complex * out = malloc(sizeof(fftw_complex) * n);
    // fftw_complex * inin = malloc(sizeof(fftw_complex) * n);
    // for(int i = 0;i < n;++i){
    //     in[i] = sin(2 * 1.3*i * 3.14 / 64) + sin(2 * 3.4 * i * 3.14 / 64);
    // }
    // write_data("inin.txt",in,n);
    // fft_r2c(in, out, n);
    // ifft(out,in,n);
    // write_data("in.txt",in,n);
    // write_data_c("out_c.txt",out,n);

    // free(in);
    // fftw_free(out);

    fftw_complex * data = malloc(sizeof(fftw_complex) * n);
    for(int i = 0;i < n;++i){
        data[i] = i+1 + (i+2)*I;
    }

    shift_freq(data,n);

    for(int i = 0;i < n;++i){
        printf("%f,%f ",creal(data[i]),cimag(data[i]));
    }
    printf("\n");

    return 0;
}